A Classical Picture of the Role of Vacancies and Interstitials in Helium-4 
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Motivated by experimental hints for supersolidity in Helium-4, we perform Monte Carlo simu- 
lations of vacancies and interstitials in a classical two-and three-dimensional Lennard- Jones solid. 
We confirm a strong binding energy of vacancies which is of the order of Lennard- Jones attraction. 
This is reminiscent of what has been found for vacancies in Quantum Monte Carlo simulations. In 
addition, we find a strong attraction and large binding energy of interstitials in two-dimensional 
simulations. This is mainly due to the formation of a pair of dislocations by clustering interstitials, 
in which minimizes the elastic deformation energy. We interpret the results in light of the properties 
of Helium-4. 



I. INTRODUCTION 



A supersolid is a conjectured phase of matter which 
possesses seemingly contradicting properties of having 
spatial long range crystalline (positional and orienta- 
tional) order as well as being superfluid at the same time. 
Potential evidence for such a supersolid phase in Helium- 
4 has been first detected as nonclassical rotational mo- 
ment of rotational inertia (NCRI) in torsional oscillator 
experiments by Kim and Chan (KC) in 2004,^'^ more 
than thirty years after the first proposals by Andreev and 
Lifshitz^ and Chester.* Based on their idea of a super- 
solid phase in a bosonic quantum crystal with vacancies 
(or interstitials), a number of theoretical explanations 
have been proposed. 



Although NCRI for Helium-4 has also been observed 
by other experimental groups, it remains contro- 
versial whether if this is a supersolid phase or if the 
interpreted supersolid phase is simply a bulk equilib- 
rium phenomenon. ^^^'^ Computer simulations have shown 
that the density matrix of Helium-4 crystal decays 
exponentially^^' and vacancies in solid Helium-4 are 
gapped^^. Non-equilibrium vacancies present in solid 
Helinm-4 attract each other and phase separate, thus 
purging the Helium crystal of vacancies. This shows 
that the ground state of a single crystal Helium-4 is not 
a supersolid^*. It is remarkable that the solid with the 
strongest quantum properties known in nature behaves 
so classically. 



Therefore, the understanding of the effective vacancy- 
vacancy interaction and interstitial-interstital interaction 
is crucial to our understanding of this phenomenon. Here 
we set to carry out a classical Monte Carlo simulation in 
order to compare with the properties of Hclium-4. We 
find that a classical modeling of Helium yields vacancy 
and interstitial binding energy that are larger than that 
in quantum case. Nonetheless, the phase remains essen- 
tially the same. 



II. SIMULATING A CLASSICAL 
LENNARD-JONES CRYSTAL 

For our classical Monte Carlo simulations we model 
the interaction between two atoms by a Lennard- Jones 
potential, where we fix units by choosing parameters suit- 
able for Helium-4:^^ 
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(1) 



with a = 2.56A and e/fc = 10.2K. 

We perform classical Monte Carlo simulation in the 
NPT ensemble of constant particle number, pressure, and 
temperature using two types of moves. The first type 
is a local particle displacement: we propose to move a 
randomly chosen particle by a random distance in the 
interval [0, a) for the simulation of vacancies in a crystal, 
and by a random distance in the interval [0, cr/4) for the 
simulation of a crystal with interstitials. The second type 
is a volume change, where we propose to vary the volume 
of the system by adjusting the length, width and height 
of the system individually within the range of ±2.5%. 

Initially, we start from a perfect triangular lattice for 
our two-dimensional (2D) or from perfect hexagonally 
closed packed (hep) lattice for our three-dimensional 
(3D) simulations. For a solid with vacancies, we ran- 
domly remove particles from the system, while for a solid 
with interstitials, we randomly add particles to the sys- 
tem. In our simulations, we define one sweep as lOOOA^^ 
trial particle moves (where N is the number of particles 
in the system) , followed by one volume-changing move in 
each dimension. While quantum mechanically Helium-4 
forms a crystal only at a pressure of above 25.46 bar, 
classically a crystal is already formed at zero pressure. 
Therefore, we use a pressure of only P = 1 bar in our 
classical simulations. In fact, the exact value of pressure 
does not matter because we are only interested in qualita- 
tive features of a classical crystal rather than quantitative 
results for Helium-4. 

We probe for attractive interactions and binding of va- 
cancies and interstitials by calculating the binding ener- 
gies and inverse compressibility. The inverse compress- 
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FIG. 1: Snapshot of a typical configuration, showing a stable 
vacancy cluster, A^o ~ 400 with three vacancies, T = O.OIK. 
Periodic boundary conditions are used. 
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FIG. 2: Binding energy Eb and inverse compressibility 
of two and three vacancies in two dimensions at T = O.OIK 
as a function of the number of particles in a commensurate 
solid A'^o- The arrows indicate on which axis to read off the 
values of the data points. 



ibility can be measured as 



(2) 



where N is the number of particles in the system. Ther- 
modynamically, phase separation will occur when k^^ < 
0. Numerically, the inverse compressibility is computed 
as 

= E{N + 1) + E{N- 1)^2 E{N). (3) 

The binding energy of a system with p vacancies can 
be expressed as 

EBiNo,p)=EiNo-p) + {p-l)E{No)-pEiNn-l) , 

(4) 

and that of p interstitials as 

Eb{No,p) = E{No+p) + {p- l)E{No)-pE{No + 1) 

(5) 

where Nq is the number of particles of the commensurate 
crystal. 



III. VACANCIES IN A CLASSICAL CRYSTAL 

For the simulation of vacancies we consider systems 
with 1, 2 and 3 vacancies respectively. To find the ground 
state, we perform an annealing process. Due to the sim- 
ple structure of energy space of vacancy systems, a rapid 
annealing schedule can be used. We start from a tem- 
perature of IK and reduce it by a factor of 2 every 100 
sweeps until we reach a final temperature of O.OIK. 

As shown in Fig. 1, vacancies have a tendency to bind 
and form a vacancy-cluster. The three-vacancy cluster 
shown to be the ground state of the classical system 
is also the most likely configuration in Quantum Monte 
Carlo simulations by Ref. 17. 
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FIG. 3: Binding energy Eb and inverse compressibility k.~^ 
of two and three vacancies in three dimensions at T = 0.017^ 
as a function of the number of particles in a commensurate 
solid iVo. 



In Figs. 2 and 3, we show two- vacancy binding energy, 
three- vacancy binding energy, as well as two- vacancy in- 
verse compressibility at temperature T — O.OIK and 
pressure P = 1 bar for two and three dimensions. 

From the Lennard- Jones potential, the classical min- 
imum interaction energy between two Helium atoms 
at zero temperature and zero pressure is exactly e = 
— 10.2K. (This is slightly higher at finite pressure or tem- 
perature.) Comparing Fig. 4a) and b), we realize that 
if two vacancies cluster, the system will be stabilized by 
the energy of one pair of atoms, which costs about -lOK 
at finite temperature and pressure. This argument shows 
that our simulation results for the binding energy of two 
vacancies in Fig. 2 and 3 are sensible. 

For a three-vacancy cluster (Fig. 4c) there are three 
pairs of interaction and the binding energy is therefore 
approximately three times as much as that of a two- 
vacancy cluster. Therefore, the value of -30K is once 
again expected. The clustering of vacancies can also be 
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FIG. 5: A typical configuration of tfiree interstitials, A'o = 256 
witli tliree interstitials, T — 10~^K. Open circles represent the 
position of the He-4 particles; dots represent an imaginary 
perfect lattice; and solid circles represent the interstitials. To 
select the interstitials, we first superpose the actual configura- 
tion onto a virtual perfect consumerate lattice. After which, 
we pair up each lattice point with the closest actual particle 
in the configuration. Finally, we remove all these pairs and 
what remains will defined as the interstitials. The final as- 
pect ratio of the configuration is 0.889 (c.f. 0.867 of a NVT 
configuration). Here, we want to emphasize the importance 
of not fixing the aspect ratio, because one will never obtain 
the lowest energy configuration once this aspect ratio is fixed. 



FIG. 4: The binding of vacancies: a) shows the configuration 
before and b) after the formation of a two- vacancy cluster. By 
binding vacancies the number of nearest neighbor atoms is in- 
creased by one and the bound state stablized. c) The three va- 
cancies form three pairs, indicated by the three dashed lines. 
As a result, a three- vacancy cluster is expected to have a bind- 
ing energy three times than that of a two- vacancy cluster. 



understood from minimizing the surface energy of a va- 
cancy cluster. The surface in contact with the crystal 
is much lower for a vacancy cluster than for individual 
vacancies. 

Similar effects have been seen in quantum Monte Carlo 
(QMC) simulations, where a binding energy of around 
one Kelvin for two vacancies and of several Kelvin for 
three vacancies has been observed^^'^**. The substantial 
reduction of the binding energy is mainly due to quan- 
tum zero point effect. Although these values are much 
reduced, they are still sufficient to cause phase separation 
of vacancies. Qualitatively, the classical picture applies 
also to the quantum Hclium-4 system. 



IV. INTERSTITIALS IN A CLASSICAL 
CRYSTAL 

Unlike the case of vacancies, where it is relatively easy 
to locate the global minimum, the interstitial system pos- 
sesses many meta-stable local minima. To search for the 
ground state, we have to perform a stimulated temper- 



ing simulation. Instead of using the much more involved 
conventional parallel tempering approach^*^"^^, we use a 
simpler method that works well in our interstitial case. 

The updates were identical to the case of vacancies, 
except that one row of atoms was fixed to prevent the 
crystal from rotating. Otherwise, the interstitials will be 
removed when the system chooses a different orientation 
by rotating. From an initial temperature of O.IK the tem- 
perature is decreased in an annealing procedure, taking 
care that the energy does not drop by more than 0.2K in 
100 sweeps. Repeated simulations (starting from differ- 
ent initial conditions) have been carried out. More than 
80% of these simulations yield the same lowest energy 
configuration, which we can safely call the global mini- 
mum. Our results again show binding of interstitials. A 
typical configuration is shown in Fig. 5. 

We observe that at low temperatures the interstitials 
will be individually inserted into neighboring rows, cre- 
ating a pair of dislocations, connected by a line of extra 
particles (see Fig. 5). This is in contrast to vacancies 
which have a tendency to cluster into compact vacancy 
clusters. 

The binding energy and inverse compressibility are 
shown in Fig. 6. Surprisingly we find even larger binding 
energies for interstitials than for vacancies. Interstitials 
cost more energy because they significantly distort the 
crystal, and much of this extra energy can be recovered 
by the binding of interstitials into a line connecting two 
dislocations as seen in Fig. 5. 
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FIG. 6: Binding energy Eb and inverse compressibility k^^ 
of of two and three interstitials in two dimensions s a function 
of the number of particles in a commensurate solid No. 



V. CONCLUSIONS 

Our classical Monte Carlo simulations of vacancies and 
interstitials in a classical solid show strong attraction and 



binding of vacancies and interstitials in a Lennard- Jones 
crystal. This is highly detrimental to the formation of a 
supersolid in the quantum mechanical case, where vacan- 
cies and interstitials are delocalized. Similar full quan- 
tum simulations of vacancies in Helium-4 show that while 
the binding is strongly reduced by quantum fluctuations, 
it is still large enough to cause vacancies to be expelled 
from the crystal instead of forming a supersolid. Al- 
though the behavior of interstitials in solid Helium has 
not been studied in detail yet (apart from a lattice model- 
ings^), it is established that they cost a finite energy gap 
A = (22.8 ± 0.7)K^'^ to create. The strong binding seen 
in the classical case combined with the larger effective 
mass for interstitials indicates that the strong attraction 
will remain in the quantum case. 

The ALPS libraries^^ were used for parallelization and 
error evaluation. The simulations were performed on the 
Hreidar cluster of ETH Zurich and the hpcpower cluster 
of HKU. We treasure useful discussions with J. Wang, G. 
H. Chen, S. Q. Shen, and P. Corboz. Financial support 
from RGC of HKSAR is also greatly acknowledged. 
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